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1. INTRODUCTION 

This paper reports on numerical experiments on the problem of moving a flexible beam. An 
optimal control problem is formulated and transcribed into a form which can be solved using semi- 
infinite optimization techniques. All experiments were carried out on a SUN 3 microcomputer. 

2. PROBLEM STATEMENT 

We consider the hollow aluminum tube depicted in figure 1. The tube is one meter long, has a 
cross-sectional radius of 1.0 cm, and a thickness of 1.6 mm. Attached to one end of the tube is a 
mass of 1 kg, and attached to the other end is a shaft connected to a motor. For simplicity, we 
assume that the torque produced by the motor can be directly controlled. Our aim is to determine 
the torque necessary to rotate the tube and bring it to rest. The maximum torque produced by the 
motor is 5 newton-meters. The equations of motion determined by application of the standard 
Euler-Bemoulli tube with Kelvin-Voigt visco-elastic damping are: 

mw„(t ,x) + Clw lxxxx (t ,x) + Elw^it ,x) - mi2 2 (t)w(t ,x) = - mu(t)x , xe[0,l] (la) 


with boundary conditions: 

w(r,0) = 0, w,(r,0) = 0, Clw^it , 1) + Elw^it , 1) = 0. (lb) 

M [0 2 (r)w(r , 1) - w„(t , 1) - «(/)] + C/Wtoxfr , 1) + Elw^it , 1) = 0, (lc) 

where w(t ,x) is the displacement of the tube from the shadow tube (which remains undeformed dur- 
ing the motion) due to bending as a function of time and distance along the tube; u(t) is the torque 
applied by the motor, and Q(r) is the resulting angular velocity (in radians per second). We shall 
denote by 6(r) the angular displacement of the rigid body (in radians). The values for the parame- 
ters in (la) - (lc) are: m = .2815 kg/m, C = 6.&9xl0 7 pascals/sec., E = 6.89xl0 9 pascals, 
/ = 1.005 x 10 r ®m 4 , M = 1.00 kg. These values are from the CRC Handbook of Material Science. 
The tube is very lightly damped (0.1 percent ). 

We consider three problems: 

P t : Minimize the time required to rotate the tube 45 degrees, from rest to rest, subject to the given 
torque constraint.. 

P 2 : Minimize the total energy required to rotate the tube 45 degrees, from rest to rest, subject to the 
given torque constraint and the maneuver time not exceeding a given bound. 

P 3 : Minimize the time required to rotate the tube 45 degrees, from rest to rest, subject to the given 
torque constraint and an upper bound on the potential energy due to deformation of the tube 
throughout the entire maneuver. 

3. MATHEMATICAL FORMULATION OF THE THREE PROBLEMS 

We will formulate the above problems P j, P 2 , and P 3 in the form of the following canonical 
optimization problem: 

P 0 : min { g°(u ,T) I g J (u ,T) < 0,j e m } 

7 e R + , u e G T 

where R+4 (yeR I7>0 ),h 4 4 { 1 , 2 , •••»»}, 

G t 4 { u e LJ0.T] I lu(r)l < 5, re [0,7] }, 


(2) 

( 3 ) 
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and g J :G t xT -> R is continuously differentiable for j e { 0, t , ,m } . We define y(u ,T) 4 

max { g'{u ,T) ) and \|/ + (u ,T) 4 max { 0,\|/(w ,T) } . 
i c a 

We shall be making use of the following functions. First, let T denote the final time. Then we 
define 

g'(u t T) 4 T (4) 

The input energy is defined as the integral of the square of the input; hence we define 
r 

g\u,T) 4 ju(t) 2 di. (5) 


Next we define 

g\u,T) 4 (0(T) - jc/4) 2 (6) 

to be the square of the angular error at the final time. We say that the tube is at rest when the total 
energy of the tube is zero. This energy is composed of the energy due to rigid body motion and 
energy due to vibration and deformation. Rigid body energy at final time is proportional to the 
square of the angular velocity. Hence we define 

g\u,T) 4 Q(T) 2 . (7) 


The kinetic energy due to vibration of the tube at time t is given by 

i 

K(t ,u) 4 y|w,(r ,x) 2 dx , 


( 8 ) 


and the potential energy due to deformation of the tube at time t is given by 

i 

P(t,u) 4 ■^■^w xx (t ,x) 2 dx . 


(9) 


We now define 

g 5 (u ,T) 4 K(T ,u), g 6 (u ,T) 4 P(T ,u). 


( 10 ) 


The tube is at rest if g\u ,T ) = g 5 (u ,T ) = g\u ,T) = 0. 

For problem P 3 , we require that the potential energy due to the tube deformation be within a 
specified range throughout the entire maneuver. This constraint has the form P(t ,u) £ f(t) for all 
re [0 ,T], where /(•) is a given positive bound function. This is a state-space constraint , and does 
not fit the canonical form Pq. However, we can replace it by an equivalent form which requires that 
we define 

T 

g 7 («,T)4 f [max{ P{t,u)-f{t) , 0)] 2 (11) 
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Then, since P(t,u) is continuous, g 7 (u ,T) = 0 if and only if P(t ,u) £ f (t) for all t e [0,7]. 

It can be shown that g J : G t xT -» R is continuously differentiable (in the topology) in u 
and / for all j e ( 1 , 2 , 7 ) . To conform with the format of problem P 0 , we relax each of the 

equality constraints by a small amount. The relaxation can be be chosen to be sufficiently small so 
as not to matter from a practical point of view. The three problems now acquire the following 
mathematical form 

P,: min { g\u ,7) I g\u ,7) - e £ 0. g\u ,T) - e S 0. g 5 (u ,7) - e <; 0, 

( 12 ) 

g 6 (u ,T) - e £ 0, ugG t } 


P 2 : minU 2 ( M ,r) I g\u,T)-T f <J0, g\u,T)-z<. 0, g\u,T)-t<. 0. 

S 5 (u ,n-e£0, g\u ,7) - e £ 0, u e G t } 


(13) 


P 3 : min [ g\u ,T) I g 3 (u ,T) - e £ 0, g*(u ,T) - e ^ 0, g 5 (u ,T)-eS 0, 

6 7 ( 14 ) 

g 6 (u,T)-e<; 0, g 7 (u ,T) - e< 0, u e G t } 

In our experiments, we set e = 10" 4 . Thus, with this relaxation, we are requiring that the final value 
of the angle 0 be in the interval [45 - 0.5 ,45 + 0.5] degrees. 

By adding an additional state variable z(t), with z(t) * 0, the above problems can be recast as 
fixed time problems on the interval [0,1] in which one has to determine not only the (time scaled) 
control u(t), but also T, the initial value of z(t) which acts as a time scale factor, and, in fact, is 
also the final time. Although the abstract form of the fixed time, scaled problems 


ueG.Te R, 


( f(u,T) I g'(u.T)Z0,jem }, 


where G 4 G x , is indistinguishable from that of the free time problem, the fixed time problem does 
not lead to the serious convergence problems that are associated with the discretization of free time 
problems. 

4. THE ALGORITHM 

To solve the above problems in fixed time form, we use an extension of the Mayne-Polak 
phase I - phase II algorithm [1], The algorithm first determines a search direction and then a step 
size to update the design parameters u( ) and T. The algorithm requires an initial guess T > 0 and 
ue G. We state this algorithm in conceptual form. 


Conceptual Algorithm 

Data: T 0 eR,u 0 eG, <xe (0, 1), pe (0, 1), p > 0 

Step 0: i = 0. 

Step 1: Compute search direction 8u, = v. -u, , 87, = x, - 7, and the optimality function 

0(u; ,7j), where v, ,x, are the solutions of the program 
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( 16 ) 


0(«i , 7, ) 4 min ( II v - u II 2 + II x - 7, II 2 + max { - py + («, , 7, ) + 
veG.xeR, 2 2 y e a 1 


xVO',.7',). 


V - u, 
x-T, 




V - «, 

x-T t 


t-VM.Ti) } } 


Step 2: Compute the stepsize X, e S » { Xe { 0, 1 ,|3,|3 2 , • • ) such that 

if \y(u, ,7, ) > 0 (at least one constraint is violated) 

X, = max ( Xe S I \| /(«, + XSu, , 7, + X87,) - y(u, , 7,) & aX0(«, ,7,) } 


(17a) 


if VO*; ,7,) £ 0, ((u,' ,7,) is feasible) 

X,- = max { Xe S I g°(Uj + XSu, , 7, + XS7,) - g 0 (“i ,7,) S aX0(u, , 7;) 
and y(u, + XSu, ,7, + X57, ) s 0 ) 


(17b) 


Step 3: Set u,' + j * u,- + X,Su,. 7, + [ = 7,' + X,87/. 

Step 4: Set i = i + 1; go to Step 1. ■ 

The function 0(-, ) is called an optimality function. It has two important properties: (i) For all 
7>0 and ueG, 0(u , 7) £ 0, and (ii) if 0(u, , 7,)<0, then (u, .7,) is not optimal and 
(v, - u, ,x, - 7,), where (v, ,X;) is the solution of (16), is a direction of descent for y if (u, ,7,) is not 
feasible and for g° otherwise. The following theorem can be deduced from the results in [1]. 

Theorem 1: If { (u, ,7, ) } is a sequence generated by the conceptual algorithm and (u, 7) is an 

accumulation point of this sequence, then 0(u, 7) = 0. ■ 

The above algorithm is called a conceptual algorithm because we cannot solve system 
(la) - (lc) exactly, and hence we cannot evaluate g J (u ,7) or Vg J (u ,7) exactly. Furthermore, since 
u is an infinite dimensional design vector, it can only be entered into a computer in discretized form. 
Hence, in practice, we must use an implementable algorithm which accepts approximations. The 
algorithm that we use adjusts integration precision adaptively, along the lines described in [2, 3 
Appendix A]. To discretize the PDE in space, we use the finite element method. Since the PDE is 
fourth order in space, it is necessary to use elements of at least second order. We have chosen Her- 
mite splines as basis elements. The input u e G is discretized in time and Newmark’s method is 
applied to evaluate the resulting system of ordinary differential equations. For a specific number of 
finite elements, p , and a number of time steps, n , the resulting discretized problem has the form: 


'n. P - min. ( gn, P (u ,T) \ gi (u ,T) 
ueG m ,Te R* 


— 


where G" 4 { u e R" I I u j I < 5 , j e a ) . 


(18) 
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The resulting problem P n p is finite dimensional and can be solved by computer. Problem P„ p 
always has a solution because the set G" is bounded. However, at first examination, it is not clear 
how solutions to P„ p relate to the solution to P 0 . Fortunately, it is possible to establish the follow- 
ing theorem which is an extension of the results in [2, 3). 

Theorem 2: Let (w„ p ,T n p ) be a solution to P„ p . If (u,7) is an accumulation point of 

{ («»./> 1 as n — > oo , p — > oo, then («, 7) is a solution to P 0 . ■ 

Implementable Algorithm 

The implepientable algorithm continues solving problem P nj) until a test indicates that both n 
and p must be incremented, i.e., the implementable algorithm increases the discretization in time and 
space adaptively. When the algorithm is far from a solution, it is less important that the partial 
differential equations be solved exactly. By using a coarse discretization in the early iterations, we 
save in two ways: the effort in solving the differential equations is smaller, and the number of 
design parameters (the size of the discretized control) is much smaller. The test for precision 
refinement monitors the progress in the reduction of \|/(m ,T), when ( u ,7) is infeasible, or in the 
reduction of g°(u ,7), when ( u ,7) is feasible. When that reduction is smaller than a parameter y> 0, 
both the number of finite elements and time steps are doubled while y is halved. The following 
theorem can be obtained by extending the results in [2, 3 Appendix A]. 

Theorem 3: Let{ (u, ,7,) } be the sequence produced by the implementable algorithm with the 

refinement criterion above. Then the discretization becomes infinitely refined as i and any 

accumulation point of { (u, ,7,) J , (u, 7), satisfies the optimality condition 0(u, 7) = 0. ■ 

5. COMPUTATIONAL RESULTS 

The results presented here are for the case in which the n 2 (r) terms are neglected in equation 
(la) - (lc). Similar results have been obtained by performing experiments for the case in which the 
Q 2 (r) terms are included. 

Problem P^ 

For simplicity, we choose the zero function as initial control and 2 for an initial value for the 
maneuver time. 

Figure 2 is a graph of the control after 150 iterations. The number of time steps is 256 and the 
number of finite elements is 48. 

Figure 3a is a graph of y n<p (u ,7) as a function of the iteration number. Figure 3b shows 
y ntP (u ,7) for the first 15 iterations. The initial discretization is 32 time steps and 6 finite elements. 
The discretization is refined at iterations 67, 99, and 123. After precision refinement, algorithm 
finds a feasible value for the control and final time for the new problem P np in only a few addi- 
tional iterations. At each refinement the value of y nj , increases. This is due to improvement in the 
accuracy of the evaluation of the partial differential equation. This increase in y nj> decreases each 
time the discretization is refined and we can show that in the limit the increase is zero. 

Figure 4 is the graph of the cost as a function of iteration number. 

Figure 5 is the graph of w(t , 1), the displacement of the tip of the tube, from the shadow tube , 
as a function of time. There is a maximum displacement of the tip of about 5 mm. This is within 
the range of validity of the Euler-Bemoulli model. The tip displacement is large between 0.36 
seconds and 0.437 seconds. 

Figure 6 is a profile of the tube deformation, w{t ,x) (see figure 1), during this interval. The 
total time for the entire maneuver is 0.7886 seconds. 
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Problem P 2 : 

Formulating the slewing problem as a minimum time problem has two drawbacks. First, the 
solution to the problem is a bang-bang control (figure 2). Bang-bang controls may be undesirable 
because they may cause premature aging of the equipment. Furthermore, bang-bang controls tend to 
excite the high-frequency modes of the system. High-frequency modes are less well modeled by 
system (la) - (lc), and it is therefore best not to excite them. Second, the simple minimum time 
formulation does not take into account the amount of energy expended in performing the maneuver. 
In certain applications, the total energy available may be limited, while the total time of the slewing 
motion is less critical. Fortunately, both of the problems arising from minimum time control can be 
mitigated by reformulating the problem. We minimize the total input energy while constraining the 
final time to be less that a specified amount. 

Figure 7 is the graph of the control produced by minimizing the total input energy while con- 
straining the final time to be less than 0.800 seconds. The resulting final time is 0.800 seconds. 
This is an increase of only 1.4 percent in the final time. The control has become much smoother 
and the total energy is reduced from 19.15 to 15.72, a reduction of 18 percent. 

Figure 8 is the graph of the control for final time being 1.00 second. This is an increase of 27 
percent in time over the minimum time case, but the total energy is reduced to 7.27, a decrease of 
62 percent. 

Problem P 3 : 

In Figure 9, curve A is the graph of the potential energy of the tube as a function of time for 
the control generated in solving the minimum time problem Pj. In problem P 3 , we have the addi- 
tional requirement to keep the potential energy, which is a measure of the total tube deformation, 
below the parabola (B) for all time. 

Figure 10 shows the optimal bang-bang control for problem P v The optimal final time for this 
case is 0.8177 seconds, an increase of 3.7 percent over the solution of problem Pj. 

Figure 11 shows the potential energy curve for the optimal control (Figure 10). 
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